Band gaps in pseudopotential self-consistent GW calculations 
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For materials which are incorrectly predicted by density functional theory to be metallic, an 
iterative procedure must be adopted in order to perform GW calculations. In this paper we test 
two iterative schemes based on the quasi-particle and pseudopotential approximations for a number 
of inorganic semiconductors whose electronic structures are well known from experiment. Iterating 
just the quasi-particle energies yields a systematic, but modest overestimate of the band gaps, 
confirming conclusions drawn earlier for CaBe and YH3. Iterating the quasi-particle wave functions 
as well gives rise to an imbalance between the Hartree and Fock potentials and results in bandgaps 
in far poorer agreement with experiment. 
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Materials whose metallic or semiconducting nature has 
not been established experimentally pose a particular 
problem for first-principles electronic structure calcu- 
lations. Density functional theory (DFT) calculations 
vifhich describe them as being metallic are not conclu- 
sive because such calculations are known to systemati- 
cally underestimate band gaps. Conventional many-body 
perturbation theory (MBPT), in which the dynamically 
screened Coulomb interaction is treated perturbatively, is 
problematic because the size of the perturbation depends 
critically on whether metallic or semiconducting screen- 
ing is used. For YH3 and CaBg, this problem was solved 
by iteration to self-consistency.^'^ Though these particu- 
lar theoretical predictions were confirmed by subsequent 
experiments, recent work casts doubt on the predictive 
capability of parameter- free MBPTJ^i^^ For this reason, 
it is important to benchmark the iterative procedure by 
applying it to a series of semiconducting materials whose 
band gaps are well documented experimentally. 

Of various many-body techniques developed in recent 
years to provide corrections to DFT descriptions of the 
electronic structures of real materials, the so-called GW 
approximation^ in which the self-energy S is approxi- 
mated by the product of the one-particle Green function 
G and the dynamically screened electron-electron inter- 
action W, has been particularly successfulASiiS Because 
GW calculations are computationally very demanding, it 
is common to make a number of additional approxima- 
tions. The Green function is usually approximated by 
a quasi-particle (QP) expression, where the QP energies 
are interpreted as a band structure. A further approxi- 
mation involves using an independent particle expression 
for G in the computational steps required to evaluate W 
and S. Generally the energy levels and orbitals obtained 
in a DFT calculation are used in this expression, most of- 
ten at the local density approximation (LDA) leveliii*i^ If 
additionally the QP wave functions are approximated by 
the LDA orbitals, then the QP energies simply become 
first-order perturbation corrections to the LDA energy 
levels. A GW calculation that involves all these approx- 
imations is referred to as "single shot" . The single shot 



GW approach yields good QP band structures and rea- 
sonable band gaps for a wide rangeSiSiiS of semiconduc- 
tors and insulators provided the reference independent 
particle spectrum, i.e., the LDA spectrum, has a band 
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gap^^^^— 

For materials such as Ge or InAs which LDA wrongly 
predicts to be nearly metallic or metallic, the single shot 
GW procedure is no longer appropriate. Iteration was 
suggested as a solution to this problem^ but in prac- 
tice the computational expense was prohibitive and in- 
stead use was made of the experimentally known dielec- 
tric screening to study such materials^i thereby introduc- 
ing an element of empiricism. For materials whose exper- 
imental character is unclear but which are described by 
the LDA as being metallic, there is no other option than 
to iterate the GW solution. Though this iterative GW 
technique has been applied successfully to YH3 with 24 
atoms in its primitive unit cell and to CaBg, its general 
validity has yet to be demonstrated. We will apply it to a 
set of materials with band gaps that range from 0.4 to 5.5 
eV: InAs, Ge, Si, GaAs, AlP, GaN and C. Note, in par- 
ticular, that the small band gap semiconductor InAs has 
a metallic LDA spectrum and that Ge has a vanishingly 
small LDA band gap. 

Computational details. Our main interest lies in de- 
veloping a practical scheme to handle complex materi- 
als, such as YH3, with large unit cells, and we have 
not attempted to go beyond the QP approximation^iiS 
For the same reason we use pseudopotentials (PP) to 
represent the ion cores, since for complex materials all- 
electron GW calculations are computationally even more 
demanding^iSiSiiS We use the real-space imaginary-time 
GW techniqueSSiSi as implemented by van Gelderen et 
alm^ Our starting point is an LDA calculation with a 
standard exchange-correlation functional, ^^'^^ normcon- 
serving pseudopotentials^'^ and an expansion of the or- 
bitals in a plane wave basis set. The LDA eigenvalues 
and orbitals are used to set up an initial Green function 
G. The dynamically screened electron-electron interac- 
tion W is then calculated using the random phase ap- 
proximation (RPA) and from G and W the self-energy S 
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is constructed. Within the QP approximation the funda- 
mental equation is 



^/.„k(r) + / dr'S^(r,r')^„k(r') 



+ y"(ir'S'=(r, r';£„k)V'rik(r') = £nkV'«k(r), 

where Vext stands for the sum of all ionic pseudopoten- 
tials, Vh is the Hartree potential and is the exchange 
or Fock potential. S"^ is the correlation part of the self- 
energy, which includes the dynamic screening and is en- 
ergy dependent. e„k and '!/'nk(r) are, respectively, the 
QP energy and wave function. The e„k are complex; the 
real parts comprise the band structure and the imaginary 
parts give the inverse QP lifetimes. 

In practice 1^ is solved by expanding the QP wave 
functions in the orthonormal basis formed by the LDA 
orbitals Vnk(r) = J2n' ^^^i' ,nkipnf^{r). Representing the 
potentials and self-energy on the same basis, is trans- 
formed into a matrix eigenvalue equation. 



(1) 1 



[hk + S|^(e„k)] C„k — £nkC„k, 



(2) 
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where (hk)„'„ = {'<fi^^\ - + i-cxt + '^if + S^^lV.'^i 
and (S^)„,„ = This equation is non- 

linear in £„k because of the energy dependence of E'^. It 
is solved by standard diagonalization and root-searching 
techniques. If the QP and LDA wave functions do not 
differ significantly, so that V'nk(r) ~ V'^k'^(r) or c„/^„k = 
(5„'„ , then Q breaks up into a set of decoupled equations 
that involve the diagonal matrix elements onljii^ii^ 



In our iterative procedures the QP energies and, op- 
tionally, the QP wave functions are used to construct a 
new Green function G. We approximate G by an in- 
dependent particle Green function, using only the real 
part of the QP energies and setting the QP weights to 
unity in order to obey the sum rule for the spectral den- 
sity. The QP wave functions are then explicitly orthog- 
onalized. After constructing a new W and S, the QP 
equation is solved again. This procedure is iterated until 
there are negligible changes in the QP energies. Since 
the first iteration starts from the LDA spectrum, it can 
sometimes be convenient to open a band gap artificially 
in this spectrum. We checked that the converged band 
gap in the QP spectrum does not depend on the size of 
the initial band gap. 

For the materials treated in this paper, all iteration 
schemes to be discussed below converge to 1 meV or bet- 
ter within 9 iterations. Spin-orbit coupling is not taken 
into account explicitly in the solids, but the end results 
are corrected a posteriori when necessary. Typical val- 
ues for the parameters used in our GW calculations are 
a real space mesh of (6x6x6) points in the unit cell, 
(8x8x8) unit cells in the interaction cell and 200 (oc- 
cupied plus unoccupied) bands iS242ii2S The results are 



Iteration on: 



• eigenvalues only 
O eigenvalues and wave functions 
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FIG. 1: Comparison of calculated and experimental band gaps 
for 7 representative inorganic semiconductors using two dif- 
ferent iterative GW schemes discussed in the text. 



extrapolated to absolute convergence in the individual 
parameters, and the absolute values of these corrections 
are summed to define a computational error bar. 

Iteration on QP energies. The simplest iteration 
scheme involves QP energies only. The QP wave func- 
tions are assumed to be identical to the LDA orbitals 
and (jJJ is solved to obtain £„k- These are used to up- 
date G, W, and E'^ before Q is solved again and the 
procedure repeated until the £„k are converged. Note 
that between iterations only the correlation part of the 
self-energy needs to be updated. 

In the case of InAs care must be taken, since in the 
LDA band structure the valence and conduction bands 
overlap. In particular, the "light hole" valence band 
mixes with the lowest conduction band leading to an 
avoided crossing around F. In our GW calculations we 
allow these two bands to "de-mix" by calculating the 
off-diagonal matrix element of the self-energy in Q and 
solving a 2 x 2 version of (01 instead of 

The band gaps calculated in this way ("eigenvalues 
only") are compared with experiment in Fig. ^ and the 
numerical values are given in the third column of Tabled 
They are seen to overestimate the experimental values by 
5 to 25%. Assuming such overestimates also hold for YH3 
and CaBg whose band gaps were calculated using the 
same approximations^s^i^ we arrive at an uncertainty 
in the predicted band gaps of < 0.2 eV, which does not 
change the prediction that these compounds are semi- 
conductors. The advantage of the iteration procedure is 
that the final, self-consistent gap does not depend upon 



3 



TABLE L Band gaps of 7 representative inorganic semicon- 
ductors. The experimental values are taken from Ref.|2^, ex- 
cept when stated otherwise. The theoretical band gaps were 
calculated using the two different iterative GW schemes de- 
scribed in the text and have had a spin-orbit splitting cor- 
rection added to facilitate comparison with experiment. The 
theoretical "error-bars" are our best estimate of the uncertain- 
ties which result from using a finite cutoff, k-point sampling 
etc. All values are in eV. 





Exp. 


QP energies 


QP energies -|- functions 


InAs 


0.4f 


0.41 ± 0.03 


1.06 ± 0.10 


Ge 


0.74 


0.80 ± 0.03 


1.72 ± 0.11 


Si 


1.17 


1.44 ± 0.04 


2.23 ± 0.09 


GaAs 


1.52 


1.62 ± 0.05 


2.66 ± 0.15 


AlP 


2.50 


2.95 ± 0.03 


3.43 ± 0.12 


GaN 


3.30" 


4.10 ± 0.05 


4.77 ± 0.31 


C 


5.49 


6.22 ± 0.04 


6.94 ± 0.17 


" Reference 



whether the initial LDA band structure corresponds to a 
metal or a semiconductor. 

Iteration on QP energies and wave functions. So far 
we have used "self-consistency" to refer to the QP ener- 
gies only, and have approximated the QP wave functions 
by LDA orbitals. While it is generally believed that the 
QP wave functions "^nk are virtually identical to the LDA 
orbitals, with overlaps |(V'^k^lV'nk)| of 99.9%fi^ there 
have been some suggestions that the differences, though 
small, might be important. By solving the full matrix 
equation Q in a single iteration, we confirm previous 
results for large overlap values but find that they only 
hold for points and lines of high symmetry. For points 
of lower symmetry, substantial deviations from unity are 
found (Fig. |2) . When the QP wave functions as well as 
the QP energies are updated the deviation becomes even 
larger. The reason why the overlaps are so large at points 
of high symmetry is quite simple: the higher the symme- 
try of a k-point, the greater the number of irreducible 
representations, resulting in a larger average separation 
of levels which can be coupled by the self-energy. Because 
the overlaps are larger than expected, we also investigate 
the effects of iterating the QP wave functions in a self- 
consistent procedure involving the solution of the full2& 
matrix equation 



The band gaps calculated by iterating the QP energies 
and wave functions are given in the fourth column of Ta- 
ble^ and are also shown in Fig.Q] They overestimatei^ 
the experimental values by between 1.0 and 1.3 eV. It- 
erating the wave functions thus has a very significant 
effect on the energy spectra, even though the final QP 
wave functions still have a large overlap (mean overlap 
of more than 90%) with the LDA wave functions. Rela- 
tively small changes in the wave functions result in large 
changes in the self-energy and in the band gaps. 
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FIG. 2: (Color online) Overlap of the QP wavefunction with 
its LDA counterpart, | ('i/'J;!k^|V'nk) | , plotted as a function of 
1/Ae, the inverse separation to the LDA state with the next 
largest overlap, for the 4 lowest QP states of silicon and 39 
irreducible k-points. At the symmetry points F, X and L the 
four occupied states are shown as shaded/coloured symbols 
and, in addition, the four lowest unoccupied states are shown 
as open symbols. 



It turns out that the Hartree and the Fock potentials, 

^^,l^«k(r')P 

respectively, 



occ „ 



r — r' 



(4) 



^^r,r') = -e'J2 

nk 



V'„k(r)Ck(r') 



r — r' 



(5) 



are especially sensitive to the wave functions. Both po- 
tentials are large, so small changes in the wave functions 
can have a large absolute effect. This is confirmed by 
fixing Vh at the LDA level while iterating and E"^ 
using the QP energies and (orthogonalized) wave func- 
tions. This leads to a band gap of 4.35 eV for silicon ! 
The Fock potential is attractive and during iteration 
its effect on the occupied states increases, whereas the 
effect on the unoccupied states decreases, which enlarges 
the band gap. This effect is partially cancelled by iter- 
ating the repulsive Hartree potential Vh, which behaves 
similarly to upon iteration, but with an opposite ef- 
fect on the band gap due to its sign. The cancellation 
however is incomplete resulting in band gaps that are 
considerably too large as seen in Table in and Fig. [2 

It is very tempting to speculate that this cancellation 
would be more efficient in an all-electron (AE) calcula- 
tion. A number of single-shot AE calculations yield val- 
ues for the band gap of silicon which are too small com- 
pared to experiment (^i^ see the third column of TablelTll 
In an even more recent paper reporting on the iteration 
of an all-electron GW scheme to self-consistency, a value 
for the bandgap of Si of 1.14 eV was foundi^S 
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TABLE 11: Comparison of the results of a number of pseu- 
dopotential (PP) and all-electron (AE) GW calculations for 
the fundamental band gap of Si, based upon the Dyson equa- 
tion (DE), or the quasi-particle approximation (QP). All val- 
ues are in eV. 



Reference Method 'Sgw{Glda) J^gw{Ggw) 



m 


PP/DE 26 lowest states 


1.34 


1.91 


\£ 


AE/DE; S diagonal 


0.85 


1.03 


M 


AE/QP 


0.90 




m 


AE/QP 


0.84 


1.14 




AE/QP 


0.92 




Present 


PP/QP 


1.17 


2.23 


Present 


PP/QP; 26 lowest states 


1.04 


1.94 


m 


Experiment 


1.17 





Because AE GW calculations are much more expen- 
sive than PP calculations, it is important to understand 
the source of the discrepancy between the two in order 
to develop computationally less demanding approximate 
schemes. If we accept that an all-electron, iterative GW 
scheme yields results in agreement with experiment then 
we must conclude that the good agreement with exper- 
iment found by single-shot plane-wave pseudopotential 
(PP) calculations is accidental and comes about from 



a cancellation of two errors. Since the single-shot AE 
calculations'*'^ *^'^^ for Si find band gaps ranging from 0.84 
to 0.92 eV, the PP scheme would appear to introduce 
an error of 0.25-0.33 eVi^fl Because the first-principles 
pseudopotentials are constructed so as to yield the same 
Kohn-Sham eigenvalues spectrum as the AE calculation, 
the error must be introduced by the pseudowavefunc- 
tions in the polarization function where it enters via the 
dipole matrix elements; this could be confirmed by ex- 
plicit comparison of AE and PP calculations. The error 
in the single-shot AE schemes are apparently almost ex- 
actly compensated on iteration^iSS In the PP schemes, 
the combination of iteration and pseudo-wavefunctions 
leads to an amplification of the errori^ 

Conclusions. We have tested two iterative, self- 
consistent GW schemes within the quasi-particle approx- 
imation, treating only the valence electrons and using 
pseudopotentials, on a range of semiconductors and insu- 
lators. Iterating the QP energies to self-consistency leads 
to an overestimation of the band gaps by 5 to 25%, which 
is sufficiently accurate for this scheme to have predictive 
power and leaves unchanged the conclusions drawn in 
Refs. l,.2j. If the QP wave functions are allowed to 
change in the iterations, the results get worse and the 
band gaps are overestimated by up to 1.3 eV. More work 
needs to be done to understand in detail the origin of this 
error. 
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This would appear to be confirmed by PP calculations in 
which the highest atomic (semi) core states were included as 
valence states in a GW calculation. At least in single shot 
calculations, inclusion of these semicore states leads to a 
decrease in the calculated band gaps.^^ However because 
they have not iterated their calculation, Tiago et al. draw 
different conclusions from this result. 



We can reconcile the value of 1.91 found by Schone and 
Eguiluz f\M with our own value of 2.23 eV by iterating 
a PP GW calculation including only 26 states; we then 
find agreement with these authors within the error bar of 
the calculation. We have no explanation for their finding a 
single-shot Si band gap of 1.34 eV. Our own value of 1.17 
(1.04 with 26 states) eV is compatible with other single- 
shot PP-GW studies. 



